require(ggplot2)
require(xtable)
require(reshape2)


delta <- rho <- seq(0.05, 1, length.out = length.out)

sim_res <- list()
for (i in 1:length(delta)) {
  fn <- paste("./nochange/corr07/res_coef/sim_coef_mse_corr07_", i, ".rds", sep = '')
  sim_res[[i]] <-  readRDS(fn)
}

## data prep for plot
mse_lasso  <- do.call("cbind", lapply(sim_res, function(x) x[,1]))
mse_bridge <- do.call("cbind", lapply(sim_res, function(x) x[,2]))
mse_ridge  <- do.call("cbind", lapply(sim_res, function(x) x[,3]))
mse_en     <- do.call("cbind", lapply(sim_res, function(x) x[,4]))

mse_lasso_dat  <- data.frame(mse_lasso)
mse_bridge_dat <- data.frame(mse_bridge)
mse_ridge_dat  <- data.frame(mse_ridge)
mse_en_dat     <- data.frame(mse_en)

colnames(mse_lasso_dat) <- colnames(mse_bridge_dat) <- colnames(mse_ridge_dat) <- colnames(mse_en_dat) <- paste("delta", delta, sep = '')
mse_lasso_dat$rho  <- rho
mse_bridge_dat$rho <- rho
mse_ridge_dat$rho  <- rho
mse_en_dat$rho     <- rho

mse_lasso_dat$method  <- "Lasso"
mse_bridge_dat$method <- "Bayes Bridge"
mse_ridge_dat$method  <- "Ridge"
mse_en_dat$method     <- "ElasticNet"

mse_full <- rbind(mse_lasso_dat, mse_bridge_dat, mse_ridge_dat, mse_en_dat)

mse_full <- melt(mse_full, id.vars = c("rho", "method"), variable.name = "delta", value.name = 'L2')
mse_full$delta <- delta[as.numeric(mse_full$delta)]
mse_full$method <- relevel(factor(mse_full$method), 'Bayes Bridge')

# mse_full$L21 <- round(mse_full$L2 / 1e5, 3)


## summary
summary.coef <- cbind(with(mse_full, tapply(L2, method, mean)),
                       with(mse_full, tapply(L2, method, sd)))
colnames(summary.coef) <- c("Mean", "SD")
summary.table <- xtable(summary.coef)
a <- print(summary.table, booktabs = TRUE)
write(a, file = "SI/Table4_column2.tex")
